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Abstract 

A low density binary mixture of granular gases is considered witliin the Boltzmann kinetic 
theory. One component, the intruders, is taken to be dilute with respect to the other, and thermal 
segregation of the two species is described for a special exact solution to the Boltzmann equation. 
This solution has a macroscopic hydrodynamic representation with a constant temperature gradient 
and is referred to as the Fourier state. The thermal diffusion factor characterizing conditions for 
segregation is calculated without the usual restriction to Navier-Stokes hydrodynamics. Molecular 
dynamics simulations are reported for comparison with the results for this idealized Fourier state. 

PACS numbers: 45.70.Mg, 05.20.Dd 
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I. INTRODUCTION 



Segregation among species in a granular mixture is of significant practical importance, and 
its description, both qualitative and quantitative, a challenging problem |l| . Many physical 
systems and state conditions are of interest, involving a range of mechanisms responsible 
for segregation across a wide range of control parameters. Here, attention is focused on a 
simple system of gas mixtures in a temperature gradient, being the latter the only mechanism 
available to induce segregation. Furthermore, each species is assumed to be composed of 
smooth, inelastic hard spheres at low density, so the inelastic Boltzmann kinetic equation 
applies One component (the impurities) is taken to be dilute with respect to the other 
(the host). In this limit the host gas is not affected by the impurities, and is described by 
its own independent Boltzmann equation. 

The macroscopic balance equations, or hydrodynamic equations, obtained from this ki- 
netic theory description have special exact solutions: the host gas has zero flow velocity, a 
constant temperature gradient and a constant pressure [3|, |j] ; the impurity gas has zero flow 
velocity, a temperature proportional to the host temperature, and a pressure proportional to 

n 

a power of that temperature [51 . Since the heat flux of the host component can be expressed 
as proportional to the temperature gradient, i.e. by Fourier's law, this state is referred to 
here as the Fourier state. The coefficients in this hydrodynamic description are defined in 
terms of coupled, nonlinear integral equations obtained from kinetic theory. In this way, an 
exact description of a mixture in a temperature gradient is obtained as the basis for explo- 
ration of segregation as a function of the temperature gradient and differences in mechanical 
properties of the two species. No explicit limitations on the magnitude of the temperature 
gradient are assumed in this analysis and all higher degree derivatives are exactly zero. This 
complements the extensive studies of thermal segregation based on applications of kinetic 

riQ 

theory at the Navier-Stokes level [6|, [TJ. 

Due to the symmetry of the problem, all spatial variation of properties occurs along 
the temperature gradient, taken to be the x axis. The segregation of impurity particles 
relative to the host gas is described by the variation of the composition no(x)/n(x), where 
no(x) and n{x) are the impurity and host densities, respectively. In the absence of the 
temperature gradient, these densities are uniform and there is no segregation. For a finite 
temperature gradient, a convenient measure of thermal segregation is given in terms of the 
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thermal diffusion factor A defined by 

d J no{x) ^d\nT{x) 
dx n{x) dx 

For the special hydrodynamic state constructed here this dimensionless factor is independent 
of X and therefore a global property of the system. For more general states, A is a local 
function and segregation properties can vary throughout the system. Here, it can depend 
only on a dimensionless form of the temperature gradient as well as ratios of the mechanical 
properties of each species (size, mass, degree of coUisional inelasticity). It can be positive 
or negative, implying that the impurity concentration increases against or along the tem- 
perature gradient, respectively (a thermal analogue of the Brazil nut and reverse Brazil nut 
effects for gravitational segregation SMl l| ) . Of course, A should vanish for mechanically 
identical host and impurity particles. 

The expression for A depends on coefficients defining the hydrodynamic fields which in 
turn are defined in terms of solutions to nonlinear integral equations following from the 
underlying kinetic theory. The kinetic theory and associated hydrodynamics for the Fourier 
state has been described elsewhere [Sj, so only a brief summary is presented in the next 
section. An analytic solution to the integral equations can be obtained as described in j^, 
and the results are extended here to arbitrary dimension d in the Appendix. More generally, 
for arbitrary clegree of inelasticity the equations are solved by a truncated Sonine polynomial 
expansion \4, |5|. The results in both cases for A and the violation of energy equipartition 
measured by means of the parameter 7 = Tq/T, are illustrated for a wide range of degrees 
of inelasticity (granular and non-equilibrium systems violate the equipartition property of 
equilibrium states and hence the steady temperatures of the two systems are different in 
general). Next, event driven molecular dynamics simulations (MD) are described for 200 
inelastic hard disks {d = 2) in a rectangular box with thermal walls to generate a temperature 
gradient. The Fourier state spatial dependence of the hydrodynamic fields for both species is 
confirmed in the bulk of the system, away from boundary layers near the walls. Comparison 
of the MD and Fourier state results for the temperature ratio 7 as a function of dissipation, 
mass ratio, and size ratio shows good agreement only for weak to moderate dissipation, 
or small mechanical differences, but only qualitative agreement is found otherwise. This 
poor accuracy is interpreted as a limitation of the Sonine approximation for this property. 
However, MD and Fourier state results for the thermal diffusion factor A show much better 
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agreement over the entire parameter space. These results are discussed further in the last 
section. 



II. THE FOURIER STATE FOR THE HOST GAS AND THE IMPURITIES 



Consider a system of smooth inelastic hard spheres {d = 3) or disks {d = 2) of 
mass m and diameter a. Inelasticity of collisions is characterized by a constant, velocity 
independent, coefficient of normal restitution a, defined in the interval < a < 1. Then, 
when two particles, i and j, with velocities Vi and Vj collide, the velocities are instantaneously 
modified to new values given by 



l + a. 
vi —{g-(T)(T, 

l + a. ^.^ 
H —{g ■ cr)a, 



(2) 



where g = V1 — V2 and ^ is the unit vector joining the centers of the two particles at contact. 
The system is supposed to be very dilute, so that the one-particle distribution function for 
position r and velocity v at time t, f{r,v,t), obeys the inelastic Boltzmann equation 
In jsl, 0] a special solution of this equation was proposed, and the state associated with it 
was called the Fourier state. It is a time-independent distribution function with gradients 
only in one direction and having a scaling form in terms of the hydrodynamic fields. 
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(3) 



In the above expressions, T is the granular temperature and n is the number of particles 
density. Both are defined from the one-particle distribution function in the usual way, 
although with the Boltzmann constant set equal to unity. It is verified that the state defined 
by the above distribution has a uniform hydrodynamic pressure p = n{x)T{x), a heat fiux 
proportional to T(x)^/^, and it exhibits a linear temperature profile, 

dT{x 



dx 



Bpa" 



(4) 



with the dimensionless constant B given by a functional of (f) js, 4|. This functional, as well 
as the function 0, was approximately determined by solving the Boltzmann equation using 
a representation of the distribution function as a truncated Sonine polynomial expansion, 
and keeping only up to bilinear terms in the coefficients ]^\. The results were shown to be 



in good agreement with molecular dynamics simulation data, at least for weak inelasticity 
(0.9 < a < 1). 

Now suppose that M additional hard spheres or disks of mass mo and diameter (Tq are 
added to the host gas in the Fourier state. For M <^ N, the effect of these "impurity" 
particles or "intruders" on the host gas distribution function is negligible, so that the results 
described above remain valid. Moreover, the dynamics of the impurities is determined 
by their collisions with the host gas particles, while collisions between intruders can be 
neglected. When an intruder with velocity Vo collides with a gas particle of velocity Vi, the 
velocities are instantaneously changed into 

, m(l + ao) , ^ 

■"o = '"0 : [90 ■ 0-) 

m + mo 

, mo(l + ao) , f-. 
= vi + ■ [go ■ (t) <t, (5) 

m + mo 

where Qq = — vi and ao is the coefficient of normal restitution for collisions between 
an intruder and a host gas particle. It is also defined in the interval < ao < 1. In 
the tracer limit being considered, the one-particle distribution function for the additional 
particles, F{r,VQ,t), obeys the inelastic Boltzmann-Lorentz equation js], [l2|. A solution of 
this equation similar to (|3]) is considered. 
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Here no is the number density of intruders and the function $ has been chosen such that the 
average velocity vanishes and the granular temperature of the impurities, defined from the 
second velocity moment of F, is To (a;) = 7T(x). The parameter 7, measuring the deviation 
from energy equipartition, must be identified (similarly to the function $) by requiring Eq. 
([6]) to be a solution of the Boltzmann-Lorentz kinetic equation. A necessary condition for 
it is that (sf 

^^^^ = Cn{x)a'-\ (7) 

with a = (a + o"o) /2 and C being a constant. The distribution function $ and the constant 
C have also been evaluated by using a truncated Sonine representation for and $ 
The distribution function F is "normal" in the context of kinetic theory, i.e. it depends on 
position only through the local densities and the temperature of the system. 

As already mentioned in the Introduction, the main focus in this paper is on segregation, 
i.e. the demixing of the granular mixture, by thermal diffusion. The amount of segregation 



can be measured by the thermal diffusion factor A defined by (jT]). The value A = indi- 
cates that no segregation induced by the temperature gradient occurs, in the sense that the 
density ratio does not depend on position. On the other hand, when A > the impurity 
concentration increases against the temperature gradient, while for A < the impurity con- 
centration increases as the temperature increases. By means of Eqs. (jlj) and (E]) the thermal 
diffusion factor can be expressed in the present 

For arbitrary values of the restitution coefficients, the equations determining the coeffi- 
cients of the truncated Sonine expansions of the distribution functions are rather involved 
(even in the bilinear approximation used in refs. 4| and js]) and they must be solved nu- 
merically. Nevertheless, in the limit of small inelasticity for the host gas, i.e. 1 — a ^ 1, it 
is possible to obtain explicit expressions for the coefficients. The results are summarized in 
the Appendix. In Figs. [Hand [21 the asymptotic expressions for A and 7 for small inelasticity 
of the host gas, as given in the Appendix, are compared with the numerical solution of the 
kinetic equations. The comparison is carried out as a function of the coefficient of inelastic- 
ity a for several values of ao, as indicated in the figures. The constant values of the other 
parameters are ci = 2, mo = m, and oq = a. It is important to realize that all the expressions 
being compared have been obtained in the Sonine approximation and neglecting products 
of three or more of the coefficients appearing in the truncated expansion of the distribu- 
tion function. Consequently, nothing can be concluded about the validity of the truncated 
Sonine approximation from those figures. It is seen that the accuracy of the asymptotic 
expression for the temperature ratio extends up to smaller values of the coefficient a than 
the corresponding expression for the thermal diffusion factor. Actually, the latter does not 
depend on a in the considered limit. This behavior was to be expected since, as pointed 
out in the Appendix, the former has been evaluated to order 1 — a and the latter to order 
(1 — a)^/^. Similar results are obtained for inelastic hard spheres {d = 3). 

III. MOLECULAR DYNAMICS SIMULATIONS 

Event driven molecular dynamics simulations of a system of inelastic hard disks {d = 2) 
have been performed to check the accuracy of the above theoretical predictions. In all the 



6 



4.5 



3.5 



2.5 



1.5 



0.5 




0.4 



0.6 

a 



0.8 



FIG. 1: (Color on line) Temperature of the impurities Tq divided by the temperature of the host gas 
T as a function of the coefficient of normal restitution of the gas particles a, for several values of the 
restitution coefficient for collisions between the gas particles and the impurities, oq, as indicated 
in the insert (the lower the curve the smaller qq). In all cases, d = 2, niQ = m, and a = ao = W. 
The solid lines have been obtained by solving numerically the Boltzmann equation in the Sonine 
approximation, while the dashed ones are the expressions in the almost elastic host gas limit given 
in the Appendix. 

simulations to be reported, the number of particles is = 200 and they are enclosed in a 
rectangular cell of sides = 4Ly = 200cr. Moreover, only a single impurity was considered, 
i.e. M = 1. The simulations started with all the particles, including the impurity, uniformly 
distributed on a square lattice and with a Gaussian velocity distribution. 



A. The Fourier state 



To generate the Fourier state, the same procedure as in ref. j4[ was employed. To in^ 
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energy into the system and induce the temperature gradient, two thermal walls 13 
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FIG. 2: (Color on line) The dimensionless thermal diffusion factor A as a function of the coefficient 
of normal restitution of the gas particles a, for several values of the restitution coefficient for 
collisions between the gas particles and the impurities, oq, as indicated in the insert (the lower the 
curve the larger ao)- In all cases, d = 2, nio = m, and a = ao = W. The solid lines have been 
obtained by solving numerically the Boltzmann equation in the Sonine approximation, while the 
dashed ones are the expressions in the almost elastic host gas limit given in the Appendix. 

located at x = and x = were considered, while periodic boundary conditions were 
used in the y direction. The temperatures of the walls were chosen accordingly with the 
theoretical prediction for the Fourier state j^]. For the cases being reported, it was found 
that the system reached, after a transient period, a steady state with gradients only in the 
x-direction and no macroscopic flow. This is consistent with the value chosen for the aspect 
ratio L^/Ly, which is outside the region in which the transversal hydrodynamic instability 
exhibited by the state we are considering shows up [l5|, ll6| for the parameters used in the 
simulations. 

Figure [3] shows the pressure, temperature, and heat flux profiles measured in a system 
with a = 0.9. The finite system generated by MD has boundary layers near x = and 
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X = Lxi where a hydrodynamic description does not hold. Outside those layers a bulk 
region in which the pressure and the scaled heat flux are uniform while the temperature 
is linear is clearly identifled. This is the region in which the theoretical predictions are 
expected to hold. For a > 0.9 a similar behavior is observed, with the size of the bulk region 
increasing as a increases. 

Once the system has reached the Fourier state, several properties of the gas and the 
impurity have been measured. The system has been divided in 20 slices parallel to the y 
axis of the same width. Moreover, the results that will be shown have been averaged in time, 
and also over 500 trajectories of the system. The emphasis has been put on the temperature 
and density proflles of both the gas and the impurity, and on the values of the parameters 
7 and A obtained from them. In order to get a systematic information of the role played by 
the several properties of the impurity (ao, ?no, ctq), only one of these properties is chosen in 
each case to differ from the gas host particles. Then, for instance, in those cases in which 
a 7^ ao, the values mo = m and ctq = a have been employed. 

The quantity 'y = Tq/T has been computed by identifying the bulk region of the system 
in which the temperature ratio is homogeneous. An example for a = 0.99 is given in Fig. 
m The different data sets correspond to vary one of the parameters of the impurity as 
indicated in the insert. In each of the cases, all the others parameters are the same for the 
host gas particles and for the impurity. It is seen that the ratio is uniform over most of 
the system. The deviation from unity is small for the differences in mass and size ratios 
considered since the system is almost elastic in all collisions. However, when there is a 
large difference between the two coefficients of normal restitution a rather strong^iolation 
of energy equipartition shows up. This is a general property of granular mixtures |8l. ll7Hl9|. 
The values of 7 reported in the following have been obtained by averaging the temperature 
ratio in the region in which it is roughly uniform. 

The coefficient A could be measured, in principle, by using its definition in Eq. ([1]). 
Nevertheless, the measurement of temperatures from the simulation data introduces more 
uncertainties than the measurement of densities. Taking into account that in the Fourier 
state the pressure is uniform, Eq. ([1]) can be transformed into 

dlnnn 

This is the expression actually used to compute A, i.e. it is obtained from the slope of hnear 
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FIG. 3: Steady dimensionless pressure, temperature and heat flux profiles for a system of hard 
disks with a = 0.9. The mechanical properties of the impurity are the same as those of the host 
gas particles. The temperature and pressure have been scaled with some arbitrary reference value, 
Tr, actually the initial temperature of the system. The symbols are simulation data, the dashed 
straight lines in the pressure and heat profiles are guides for the eye, and the dashed line in the 
temperature profile is a linear fit of the data in the bulk of the system. 
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FIG. 4: (Color on line) Profile along the x direction of the ratio between the temperature of the 
intruder Tq and the temperature of the system T for a system with a = 0.99. 

profile of In Hq as a function of In n. An example is provided in Fig. [5l The three sets of data 
correspond to the intruder differing from the host gas particles in the mechanical property 
indicated in the insert. 

B. Comparisons for 7 and A 

In the following, results will be restricted to the case a = 0.9. Figures [6] and [7] show the 
behavior of 7 as a function of {m = rriQ and a = ctq) and of mo/m {a = and a = ao), 
respectively. A good agreement between the theoretical predications and the simulation data 
is observed, although quantitative discrepancies in 7 appear for small values of ao, i.e. very 
inelastic intruder. Figure [7] demonstrates the existence of segregation. The thermal diffusion 
coefficient is positive when ao < 0.9, i.e. the impurity concentration is higher at the colder 
part of the system. For larger values of ao, segregation occurs in the opposite direction; the 
impurity concentration is higher in the hotter part of the system. Consistently, the value 
of ao for which the direction of the segregation effect changes is also the value at which the 
temperature of the impurity equals the temperature of the host gas, since at this point the 
intruder is equivalent to the gas particles. 
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FIG. 5: (Color on line) Logarithm of the (dimensionless) intruder density as a function of the 
logarithm of the host gas (dimensionless) density, for a system with a = 0.99. 

The dependence on the mass ratio of both 7 and A is shown in Figs. IH]and|9l respectively. 
Now it is observed that while the theory accurately predicts the thermal diffusion factor, 
it clearly fails to describe the breakdown of the energy equipartition. When the mass mo 
of the intruder is larger that the mass m of the host gas particles, the simulation results 
indicate that Tq/T grows rather fast with nio/m, while the theory predicts a weak decrease 
remaining smaller than unity. Violations of energy equipartition in homogeneous granular 
mixtures much stronger than predicted by the existing kinetic theory models have been 



reported lately 



19|. 



Finally, the dependence of 7 and A on the diameter ratio ao/cr is given in Figs. [TOjandfTTl 
respectively. Again, the theoretical prediction for A can be considered as quite satisfactory, 
but not that for the temperature ratio. 



IV. DISCUSSION 



A theoretical description of thermal segregation has been explored under very controlled 
conditions. The fundamental assumptions are the validity of Boltzmann kinetic theory for 
a binary mixture of inelastic, smooth, hard spheres or disks, and a special exact steady 
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FIG. 6: Ratio between the impurity temperature To and the host gas temperature T in the bulk 
of the system, once the Fourier state has been reached, as a function of the coefficient of normal 
restitution for the collisions between the intruder and the host gas particles. The coefficient of 
normal restitution for the gas is a = 0.9. The mass and diameter of the intruder are the same as 
those of the host particles. The symbols are MD simulation results, the solid line the theoretical 
prediction from the Boltzmann equation using a truncated Sonine expansion, and the dashed line 
the weak inelasticity host gas limit given in the Appendix and discussed in the text. 

"normal" solution. The kinetic theory is limited to low density gases, which excludes many 
experimental conditions of interest. However, it provides a useful testing ground for explor- 
ing the rather large parameter space of binary mixtures. The normal solution is restricted to 
conditions for which all space and time dependence can be captured by the hydrodynamic 
fields, and hence excludes boundary layers as discussed above. The only segregation mecha- 
nism considered in the special solution here is a thermal gradient. However, and important 
feature distinguishing the analysis from others in the literature is the absence of any explicit 
restriction to small temperature gradients, i.e. a description of thermal segregation outside 
the Navier-Stokes limit. Still, the absence of solutions for extreme dissipation (see Figs. [T] 
and [21) suggests an implicit limitation of the normal solution. 

As expected, the MD simulation of a finite system with externally imposed temperature 
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FIG. 7: The dimensionless thermal diffusion factor A for the same system as in Fig. [6l 

gradient does not give the Fourier state exactly, due to boundary layers. However, in the bulk 
of the system the Fourier state is confirmed to good accuracy. The comparison of MD and 
theoretical predictions for the thermal diffusion factor A shows good agreement across the 
parameter space of ao, mo, ao. In all cases conditions are found for increased composition 
of the impurities both along and opposite the temperature gradient. The corresponding 
comparison for the temperature ratio 7 is less satisfactory, particularly for variations of 
mo/m and ao/a. It might be taken as a signature of the breakdown of the low order Sonine 
expansion truncation, but then it is puzzling why this same approximation should work so 
well for A. This discrepancy between theory and simulation for 7 remains unexplained at 
this point. 
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FIG. 8: Ratio between the impurity temperature To and the host gas temperature T in the bulk 
of the system, once the Fourier state has been reached, as a function of the ratio of the masses 
niQ/m . The coefficients of normal restitution for the gas-gas and for the gas-intruder collisions 
are a = ao = 0.9. The diameter of the intruder is the same as the diameter of the host particles. 
The symbols are molecular dynamics simulation results, the solid line the theoretical prediction 
from the Boltzmann equation using a truncated Sonine expansion, and the dashed line the weak 
inelasticity host gas limit given in the Appendix and discussed in the text. 

Appendix A: Nearly elastic host gas limit 

The distribution functions of ([3]) and (l6l) must be solutions to the Boltzmann and 
Boltzmann-Lorentz equations, as described in [5]. They are approximated by the truncated 
Sonine expansions. 
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FIG. 9: The dimensionless thermal diffusion factor A for the same system as in Fig. [8l 



with the constants appearing in these expressions determined from moments of the kinetic 
equations. For 1 — a; ^ 1 a systematic determination of these coefficients has been carried 
out for (i = 2 in the Appendix of j^. The general results for arbitrary dimension are given 
here: 

floi ~ ^01 = 0{1 - a), 

boi = 6io + C(l - a), 
"2(^(1 - a) T'/' 



^10 



Oil -a), 



B, 



01 



d-l 

5io = 5oi + 0(l-a), 
16(rf - l){a/aY-^ + 27h{2 - hf'^ 



C(l - a). 



(A3) 
(A4) 

(A5) 
(A6) 
(A7) 



/ii/2 [2A{d + 2) - ^{d + 8)/i + 27/^2] 
All the dependence on the mass ratio m/rrio and the coefficient of restitution occurs 
through the parameter 
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Moreover, the coefficients B and C characterizing the hydrodynamic profiles defined in Eqs. 
(jl]) and (171) have the expressions 
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FIG. 10: Ratio between the impurity temperature To and the host gas temperature T in the bulk 
of the system, once the Fourier state has been reached, as a function of the ratio of the diameters 
(To/a . The coefficients of normal restitution for the gas-gas and for the gas-intruder collisions are 
a = aQ = 0.9. The mass of the intruder is the same as the mass the host particles. The symbols are 
molecular dynamics simulation results, the solid line the theoretical prediction from the Boltzmann 
equation using a truncated Sonine expansion, and the dashed line the weak inelasticity host gas 
limit given in the Appendix and discussed in the text. 
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while the temperature ratio reads 
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(All) 



Notice that all the quantities above have been explicitly evaluated up to order (1 — a)^/^, 
with the exception of the temperature ratio 7 that has been computed up to order 1 — a. 
Finally, the thermal diffusion coefficient A is obtained by substituting Eqs. (1A9p and (lAlOp 
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